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Abstract. The evolution of three-dimensional, large- 
scale, magnetic fields, in a galactic disk is investigated 
numerically. The N-body simulations of galactic dynamics 
are incorporated into the kinematic calculations of induc- 
tion equations to study the influence of non-axisymmctric 
gas flows on the galactic magnetic field. The time- 
dependent gas velocity fields are introduced as input pa- 
rameters for the MHD-simulations. The effects of inter- 
stellar turbulence are given as the turbulent diffusion of 
the magnetic field. Our principal concern is to check how 
dynamical evolution of the galactic gas affects the global 
magnetic field structure and intensity. 

We have found that the magnetic field responds sen- 
sitively to changes in the gas velocity field, and even 
slight variations of the dynamical parameters, such as the 
gas mass/total mass ratio results in nonuniform intensity 
structures, i.e. magnetic arms. The magnetic lines of force 
are well aligned with spiral arms and bar due to compres- 
sional flows in these features. In the inter-arm regions the 
areas with magnetic vectors going opposite to the main 
magnetic spirals are present. Moreover we could identify 
Parker-like intensity features perpendicular to the disk 
plane. Obviously even weak perturbations of the velocity 
field already produce almost instantaneous changes in the 
local magnetic field geometry. This proves the importance 
of our approach to use the velocity fields resulting from 
global dynamical instabilities for magnetohydrodynamical 
simulations. 
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1. Introduction 

Spiral galaxies are magnetised disks. The evolution and 
origin of their large-scale magnetic field is supposed to be 
related to induction processes, via motion of highly con- 
ducting plasmas (for recent reviews see Wielebinski and 
Krause 1993; Kronberg 1994; Beck et al. 1996). This rela- 
tion between the plasma motion and the magnetic field is 
well known as a dynamo. The conventional dynamo the- 
ory contains a large-scale velocity typically given by the 
differential rotation of the disk plus some turbulent con- 
tribution which describes the turbulent diffusion and the 
helicity. Since the helicity is of the order of 1 km s~^, it 
is not an observable quantity and its value and physical 
origin can only be determined from theoretical considera- 
tions (Hanasz and Lesch 1993; Beck et al. 1996). 

The dominant term for magnetic field amplification is 
the large-scale velocity field, which is a really observed 
quantity. The large-scale velocity field of a spiral galaxy is 
the product of the differentially rotating disk, an almost 
rigidly rotating galactic bulge, and non-axisymmetric ve- 
locity components, which are the result of dynamical in- 
stabilities (e.g. Binney and Tremaine 1987 and references 
therein). Whereas the energy density of the large-scale 
magnetic field is a negligible portion of the rotational en- 
ergy density of a galactic disk (VFmag (x 10^^'^ergcm^'^ 
{B - S^iG), Wrot oc 10-1° erg cm-3 (w^ot ~ 200 km s'^)), 
it is comparable with the kinetic energy density of the 
perturbations driven by dynamical instabilities like spi- 
ral arms and bars (Wdyn ~ lO-^^ergcm"^ [vdyn oc 
lOkms-i)). Since the dynamical instabilities present 
global responses of the self-gravitating, dissipative compo- 
nent of the interstellar medium in disk galaxies, they will 
have a considerable influence on the structure and evolu- 
tion of large-scale magnetic fields in disk galaxies (Chiba 
and Lesch 1994). 
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A tight relation between the galactic velocity struc- 
ture and the _B-field is indeed clearly visible in many ra- 
dio observations: recent high frequency (Faraday rotation 
negligible) observations of radio polarisation from nearby 
galaxies: M83 (Neininger et al. 1991), M51 (Neininger 
1992, Neininger and Horellou 1996) and NGC6946 (Ehle 
& Beck 1993) show that the large-scale galactic mag- 
netic field is aligned with the spiral arms and bars, which 
suggests the presence of the dynamical influence of non- 
axisymmetric gas flows on the galactic magnetic field. 
High angular resolution interferometer observations sup- 
port this hypothesis, showing that the galactic magnetic 
fields follow the dust lanes in M31 (Beck et al. 1989) - 
regions with a high ratio of gas compression. Narrow lines 
of highly polarised emission have been found at the inner 
edges of the spiral arms in IC 342 (Krause 1993). These 
facts could mean a close dynamical connection between 
magnetic fields and interstellar gas flows in the spiral arms 
and bars. Lower spectral resolution radio observations, 
which are more sensitive to extended regions of polarised 
emission, show the presence of such areas also in the inter- 
arm regions (Beck & Hoerncs 1996). The magnetic field 
vectors there show a spiral structure which has the same 
pitch angle as do the spiral arms in NGC 1566 (Ehlc et al. 
1996), MSI (Krause et al. 1989) and NGC 4254 (Soida et 
al 1996). 

The appearance of the dynamical instabilities does 
not depend on an existing magnetic field (Binney and 
Tremaine 1987). Thus, we have to investigate the mag- 
netic response to a large-scale non-axisymmetric instabil- 
ity. Since the dynamical instabilities are inherently time- 
dependent, the magnetic field will also vary with time. We 
will present time-dependent 3D dynamical simulations of 
galactic disks, containing stars and molecular-gas clouds. 
The excited instabilities will lead to a time dependent ve- 
locity field, which will be inserted into the induction equa- 
tion. 

It is the aim of this contribution to present and discuss 
the results of this procedure. We consider the following 
problems: (1) What is the structure and evolution of the 
large-scale galactic magnetic field under the influence of 
spiral and bar structures present in a galactic disk? (2) 
To what extent could the resulting magnetic field explain 
the observed spiral pattern of magnetic fields in nearby 
galaxies? 

With the same motivation Otmianowska & Chiba 
(1995, Paper I) have shown from their simulation that the 
magnetic field follows closely the interstellar gas motions, 
resulting in structures aligned with the spiral arms and the 
bar, a configuration which is in good agreement with ob- 
servation (Ehle & Beck 1993, Neininger 1992, Neininger et 
al. 1991). However, this effect is visible only near the spiral 
arms and rather than in the inter-arm regions, where the 
derived magnetic field vectors go outward from the galac- 
tic center. This is not observed in actual galaxies. This 
fact is due to certain simplifications of the first model: 



the given gravitational potential of stars and the single 
(not evolutionary) velocity field of the interstellar gas. The 
present project constitutes fully dynamical evolution of 
the gas velocity field over lO'' yr in 3D and simulates two 
components of interstellar medium: stars and molecular 
gas. The mean-field dynamo mechanism is not considered 
in our calculations. 

The outlay of the paper is as follows: Sect. 2 is de- 
voted to a description of the simulation methods. The in- 
put model parameters are presented in Sect. 3, and the 
results of the simulations appear in Sect. 4. Conclusions 
and discussion are given in Sect. 5. 

2. Models and numerical methods 

First, to obtain realistic velocity fields for calculations of 
magnetic field evolution (solving the induction equation), 
we perform a three-dimensional N-body code, using the 
particle-mesh scheme. Next the resulting velocity fields of 
the gaseous component are used to calculate the three- 
dimensional large-scale galactic magnetic field. The N- 
body simulations obtain two components: stars and molec- 
ular clouds. The clouds represent the major mass compo- 
nent of the interstellar gas. Since typically the gas mass is 
only 10% of the total mass, the gas reacts like a test par- 
ticle to the gravitational potential built up by the stellar 
and dark matter background. The importance of the gas is 
emphasized given by the fact that only in the presence of 
the dissipative gas component are dynamical instabilities 
long-lived features in a disk (Binney and Tremaine 1987). 
Without a considerable amount of gas, like in SO-galaxies, 
no spiral structure or bars appear; this is supported by 
many observations and N-body simulations (Sellwood and 
Carlberg 1984; Kennicutt 1989). The magnetic field is in- 
deed connected with the gaseous component of a galaxy 
and the instabilities responsible for structure will produce 
a non-axisymmetric velocity field, which in turn will have 
a significant influence on the magnetic field. 

2.1. N-body numerical scheme 

We use self-consistent 3D simulations for the modelling of 
the galaxy. The code is a 3D version of the 2D code of 
Combes & Gerin (1985) kindly supplied by F. Combes. 

The disk of the galaxy consists of stars (iVg ^ 38 000) 
and clouds (iVg < 19 000). The particles are interacting 
gravitationally with each other and are embedded in a 
bulge and a halo. The gravitational potential is calculated 
from the density via fast Fourier transforms (FFT) using 
a 128 X 128 x 64 grid (with a useful grid of 64 x 64 x 32). 
We calculate the density at each spatial point and at each 
integration step. 

Collisions between molecular clouds are highly inelas- 
tic, and a particle model can provide the best control of 
the corresponding energy dissipation. To simulate a colli- 
sion in the cloud ensemble, we adopt the same scheme and 
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local processes as described by Casoli & Combes (1982) 
and Combes & Gerin (1985). The cloud ensemble is rep- 
resented by its mass spectrum (masses between 5 • 10^ and 
5 • 10^ M0). We also take into account that giant molec- 
ular clouds (CMC) have a finite lifetime. After about 10^ 
years, the clouds are destroyed by tidal forces and/or star 
formation (stellar winds, supernova explosions etc.). To 
describe this destruction numerically we explicitly focus 
on the CMC (mass of one CMC Mgmc > 2 • lO^M©). 
After 4 • 10^ years, cloudlets are reinjected into inter- 
stellar space with a steep power-law spectrum and a 2- 
dimensional random-velocity distribution. A star particle 
is added at this location to the stellar disk. 

We use two analytical components to simulate the 
bulge and dark halo. These components are Plummer 
spheres $buigc(r) = -GAf^ • (r^ + 62)-i/2 and $haio(r) = 
—GMh ■ {rf + a^)^^/^. The halo mass within disk radius 
Td is M/i, a is its scale length and the bulge mass is M^ 
(within radius rh) with a scale length b. 

The velocity distribution of the particles is then chosen 
to stabilise the disk with respect to axisymmetric distur- 
bances. For a stellar disk to be stable against these distur- 
bances, the Toomre parameter Q ~ ^ has to be larger 
than 1. Here v denotes the velocity dispersion of the dif- 
ferent components (gas or stars), and S is their surface 
mass density, n is the epicyclic frequency. When Q is in 
the range 1 < Q < 4 , the stellar disk is, however, unsta- 
ble with respect to non-axisymmetric perturbations and 
in particular to bar formation (Hohl 1971). 

One can define the Toomre parameter Q (Toomre 
1964) separately for gas and stars: 



Qg 
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but these are only indicative, since there is no indepen- 
dent stability criterion. The coupling between gas and 
stars favours instabilities: both components could be un- 
stable together in situations where they would be stable 
separately (Jog & Solomon 1984; Jog 1992). The initial 
circular velocity is computed with the N-body integrator, 
in order to take into account softening and disk trunca- 
tion; then the initial rotational velocity is computed from 
the Jeans equations, taking into account the velocity dis- 
persions. 

The present code has been tested in many simulations 
(e.g. Combes et al. 1990). The collision parameters (num- 
ber of particles, elasticity of collisions, etc.) have been var- 
ied around the values chosen here, without any effect on 
the results (e.g. von Linden et al. 1996). 

2.2. Solving the induction equation 

The evolution of the large-scale magnetic field is simulated 
with the three-dimensional ZEUS code, in which the re- 
cently developed constrained transport (CT) algorithm is 



implemented for the numerical evolution of the compo- 
nents of the magnetic field in three dimensions (Stone & 
Norman 1992). We use only the part of the magnetic field 
evolution, where the flux-freezing approximation applies: 



dB/dt = rot('U x B), 



(2) 



where v is the velocity of gas and B is the magnetic in- 
duction. The numerical method fulfills automatically the 
condition divB — 0. Diffusion effects are approximated by 
"pulsed-flow" method (Bayly & Childress 1989), where 
first evolution of perfectly frozen in magnetic field is con- 
sidered. Then after a few tens of time steps the effects of 
the magnetic diffusion are introduced using the solution of 
the diffusion equation with an application of the Green's 
function (Otmianowska-Mazur et al. 1992). In practice 
this is equivalent to convolving the whole grid of mag- 
netic field components with a Gaussian function having 
dispersion ^/(2r]AF), where rj is the magnetic field diffu- 
sion coefficient. The convolution method also satisfies the 
divergence-free condition of the magnetic field. Computa- 
tions with the solution of the induction equation, including 
magnetic diffusion, are in preparation. 

To incorporate the particle-based flow into the grid- 
based code, we adopt the approximation of a spline func- 
tion. The method is described in detail in Paper I (Paper 
I, eq. no. 4 & ref. therein). 

The computations have been performed using a Cray 
Y-MP in HLRZ Jiilich and Convex Exemplar SPP-1000 
Series in ACK CYFRONET-KRAKOW. 

3. Input models parameter 

3.1. N-body model 

We simulate non-interacting disk galaxies using a 3D code 
to study the effect on the vertical extension of the disk. 

The time step of the N-body simulation is 10® yr. The 
time step for collisions has not to be so frequent, since the 
collision time-scale is of the order of 10^ yr. 

The stellar particles are initially distributed in a 
Toomre disk (Toomre 1963) of mass Mj within disk ra- 
dius r J and scale length d. The gas clouds are distributed 
in an exponential disk, with scale-length 12 kpc, and mass 
4 • 10^ M0, i.e. corresponding to 4.4 % of the disk mass. 

We made several different simulations by varying the 
input parameters for the disk and the bulge. Depend- 
ing on the input parameters, the simulations produce bar 
and/or spiral structures or only axi-symmetric structure 
like rings. We are primary interested in the qualitative re- 
sults. Thus we present here only barred galaxies and we 
focus on the development of strong bars. Friedli & Benz 
(1993, 1995) explore the influence of different initial condi- 
tions of the gas, disk, bulge, and halo mass on the forma- 
tion of bars. For this reason, we concentrate on two basic 
experiments which differ in the input parameters (sec Ta- 
ble 1). 
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Model 


I 


II 


Mass in IO^^Mq: 






- dark-halo mass Mu 


9.6 


9.6 


- disk mass Md 


7.2 


7.7 


- gas mass Mg 


0.4 


0.4 


- bulge mass Mj, 


4.8 


4.3 


- mass ratio {Md + Mg)/Mtot 


0.35 


0.37 


Scale length in kpc: 






- dark halo h 


15 


15 


- bulge 6 


1.85 


2 


- disk d 


6 


6 



ing on the grid spacing (Hawley et al., 1996). Our concern 
is to clarify how the numerical diffusion affects our mag- 
netic simulations. To solve this problem, the calculations 
for model I are performed for two cases (see Table 2): 51 
grid points in the X and Y directions with step size of 600 
pc (la) and 101 grid points with spacing 300 pc (lb). Both 
experiments are made with no additional diffusion (77 — 0) 
to check how the grid spacing effects the truncation. 



Table 1. Input parameters for model I and IT 



3.2. Magnetic field evolution model 

The model parameters adopted in the computations are 
summarised in Table 2. Our experiments are divided into 
two main groups: magnetic field evolution in a galaxy pos- 
sessing gravitationally strong (I) and weak (II) bar distur- 
bance (Table 2). 

The kinematic evolution of the magnetic field is solved 
in a 3D rectangular grid of points, where the XY plane 
is the galactic plane and the Z-axis is the axis of galactic 
rotation pointing in the northern direction. The rectan- 
gular size in the X and Y directions is 30 kpc with the 
grid interval of 600 pc (51 grid points) or 300 pc (101 
grid points), depending on the experiment. The height of 
the model galaxy is 1 kpc with spatial resolution of 50 
pc (21 grid points) (Table 2). The adopted time step is 
2 • 10^ yr. Computations with a smaller time step, 10** yr, 
give identical results. The input magnetic field configura- 
tion is purely axisymmetric, going in the anti-clockwise 
direction, although a more general structure would be the 
superposition of several magnetic modes (see Paper I for 
discussion). The initial field is toroidal, its intensity de- 
creases toward the galaxy center (from R = 4 kpc) to 
avoid numerical instabilities and at radii larger than 4 
kpc is taken as a constant. At i? > 12 kpc it decays expo- 
nentially. The z-component is initially zero. A calculation 
with more realistic seed fields is in preparation. The inten- 
sity of the initial field is 10~^ G, since fields of this order 
are widely observed in spiral galaxies and the actual value 
of the initial field does not change the results (Paper 1). 
We adopt so called "outflow" boundary conditions (Stone 
& Norman, 1992). Physically this means that magnetic 
flux outflow from the galaxy is present. Such boundary 
conditions result in a strong diminishing of the magnetic 
energy density during all the experiments. The decrease of 
magnetic energy is connected with the effective magnetic 
Reynolds number in the ZEUS code. The magnetic field 
decay simulations made by Hawley et al. (1996) show that 
the diffusivity provided only by numerical truncation er- 
rors gives Reynolds numbers of order lO"' and 10**, depcnd- 



Fig. 1. a: Model I at time step 8 • 10 yr. The upper picture 
shows the total gas (right) and 9500 stellar (left) particles of 
the X — y plane. The lower picture shows the belonging y — z 
plane of the simulation. 



We adopt a diffusion coefficient 77 based on the turbu- 
lent motions of the interstellar gas. The value of 77 changes 
from 10^^ cm^/s (models Ic and lib) by 3.6 • 10^^ cm^/s 
(le) to 10^^ cm^/s (Id) and is two or three-dimensional 
(see the discussion below). The computations can last no 
longer than about 6 • 10* yr, because of the significant 
decrease of magnetic energy due to numerical dissipation 
and the outflow boundary conditions. In order to check 
late as well as early epochs of the velocity field evolution 
the calculation are performed from 7 ■ 10* to about 10^ yr 
for model I and from 2 • 10* to about 8 • 10* yr for model 
II. The computations have also been done for two dimen- 
sional diffusion working only in the XY plane (Table 2, Id 
and le) to study the possibility of much smaller magnetic 
dissipation in the Z direction than in the galactic plane. 
The diffusion coefficient value in this case is 10^^ cm^/s 
(Id) or 3.6 • 10^6 cmVs (le). 

The choice of two diffusive models needs some more 
explanation. In our philosophy, the magnetic field evolu- 
tion is a direct consequence of the dynamical instabilities 
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Fig. 1. b: Model I at 8.5 • 10* yr 



Fig. 1. d: Model I at 10^ yr 



Model dx{= dy) dz 
pc pc 



iVx 



iVy-A^Z 



cm^/s 



la 


600 


lb 


300 


Ic 


600 


Id 


600 


le 


300 



50 51 -51 -21 

50 101 -101 -21 

50 51 • 51 • 21 10^** 

50 51 • 51 • 21 lO^^planar 

50 101 • 101 ■ 21 3.6 ■ lO^'^planar 



Ila 


600 


50 


51 -51 -21 





lib 


600 


50 


51 -51 -21 


10"« 



Table 2. The magnetic field evolution model parameters 



scale Ac and velocity dispersion Vc can be defined with the 
critical value Q = \ (Lesch 1993) 



A. 



27r2G'S 



and 



ttGS 



(3) 



Fig. 1. c: Model I at 9 • 10* yr 



driven by non-axisymmetric perturbations of the gravita- 
tional potential. The molecular clouds are very sensitive 
to such perturbations. The formation of the clouds also 
needs an instability, which produces the local gravitational 
collapse. In disk systems, such a local process can be de- 
scribed by a global stability criterion, depending on the 
disk rotation, i.e. the epicyclic frequency and the surface 
mass density E (Binney and Tremaine 1987). In terms of 
the Q- value (see Eq. 1) the condition for axisymmetric 
collapse can be rewritten as Q < 1 and a typical collapse 



As an initial condition we assume that the disk already 
contains molecular clouds. That is, the disk is kept on 
the border of gravitational instability, creating a cloudy 
medium with both cloud sizes and separation of order 
Ac. These clouds should provide the angular momentum 
transport by a two dimensional random walk with a step 
length Ac and velocity Wc, i.e. the viscosity is (Lin and 
Pringle 1987) 



V - AcWc -^ >^1^{R)- 



(4) 



We note that from the onset condition Q < 1, one can 
define a critical surface mass density Sc. If the theory is 
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Fig. 2. a: Simulation for model II at time step 2 ■ 
caption Fig. la). 



10** yr (see Fig. 2. b: Model II at 3 • 10* yr 



correct then for S > Sc, cloud formation and star forma- 
tion should start. Indeed, Kennicutt (1989) and Elmegreen 
et al. (1994) found that in disk galaxies of all Hubble 
types, the size of star forming regions is of the order of 
Ac, and that the star formation starts at a radius where 
E ~ Ec. They concluded that the star forming threshold 
is associated with the onset of the large-scale gravitational 
instability in the galactic gas disks. 

Thus, the gravitational disk-instability may serve 
as the ultimate pool for turbulence in the interstellar 
medium. We take the gravitational viscosity equivalent 
to the magnetic diffusivity. According to Eq. 4 a typical 
galactic disk model gives a characteristic value of 



10^^ cm^s-i 



A 



lOOpc 



n 



10 



-15, 



(5) 



(Lesch 1993). Since the cloud formation is restricted to the 
disk, the diffusivity acts only radially in the disk plane, 
provided the disk activity is low, so that no considerable 
gas flow into the halo is expected. This is obviously the 
case in most disk galaxies, which show no extended radio 
halos (Dumke et al. 1995). 

4. Results 

4.1. Model I 

During the first time steps the disk exhibits a strong spiral 
structure. 

Since the stars are collision-free the stellar spiral arms 
are much broader than the gaseous spirals. The winding 
of spiral arms can lead to the formation of pseudo-ring 



Fig. 2. c: Model II at 5 • 10* yr 



structures (Fricdli & Bcnz 1993) which is seen in the very 
outer part of the disk. 

After 3 • 10* yr a 3-armed spiral structure appears. 
The mass transport, initiated by the spiral arms, leads to 
a strong mass concentration in the central region, which 
changes the behaviour of the disk. First, the mode changes 
into a 2 armed spiral and few time steps later the mass 
concentration allows the formation of a very massive and 
long bar. This bar reaches its maximal strength within a 
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Fig. 3. a: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model lb 
with no magnetic diffusion at t = 8 • 10 yr (top) and at t = 8.5 • 10* yr (bottom) 



few time steps. After 8 • 10® yr the bar is 25 kpc long and 
ended at the bar-corotation radius (12.5 kpc) (see Fig. la). 

The gas bar is very thin in comparison to the stel- 
lar one (see also the simulation CSGl from Combes & 
Elmegreen 1993). During the bar evolution a 2 armed spi- 
ral mode is excited in the outer disk radii is starting at 
the bar-corotation radius. 

The stars are collision-free so their Q-parameter is 
high. Therefore the stellar disk is thicker (dynamically 
hotter) than the gas disk, which can be seen in the y — z di- 
agram (see Fig. 1). The evolution of the stellar disk is very 
similar to the simulation made by Gcrin et al. 1990. We 
found the same peanuts shape bulge in the y — z diagram, 
which is formed due to the bar formation. 



The bar pattern speed decreases slowly (Fig. lb and 
c), but the bar is visible until the end of our simulation 
(1.4 • 10^ yr) (Fig. Id). We present here the same time 
steps of the dynamic evolution as we show for magnetic 
field experiment (see Sect. 4.3). 



4.2. Model II 

Initially this model is similar to the case I. A spiral struc- 
ture appears earlier in the disk than in model I, (after 
2-10® yr) (see Fig. 2a). the arms are more open than in 
the first simulation (Fig. 2b) . The mass concentration in 
the inner disk is higher than in model I. For example after 
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Fig. 3. b: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model lb 
with no magnetic diffusion at t = 9 • 10 yr (top) and at t = 9.5 • 10* yr (bottom) 



9 • 10^ yr model II has accumulated 30% more mass in the 
inner 1 kpc than model I. 

After 4 • 10* yr a bar is formed in the disk. This bar 
increases its length during the next 10 time steps (see 
Fig. 2c). The bar is less massive and shorter than in the 
first simulation and has his longest extend (20 kpc) after 
6 • 10* yr. In the outer part of the disk a two armed spiral 
pattern is situated like in the simulation made by Gerin 
et al. 1990 . This has been verified by an analysis of the 
potential via FFT like Sparke & Sellwood (1985). The 
ring-like structure at the bar belongs to orbit families at 
Lagrange point (Contopoulous & Papayannopoulos 1980). 

In the z-direction one can see again that the gas disk 
in the inner part gets thinner during this evolution. In the 



outer part of the disk r>30 kpc the disk gets thicker (3 
kpc, Fig. 2 y-z diagram). 

During the next time steps the bar shrinks and more 
mass is concentrated at the bar. We follow the bar evolu- 
tion for 1 ■ 10^ yr. 



4-3. The 3D evolution of magnetic field: model I 

The calculations with the strong bar disturbance (model 
I) have been done for five sub-cases: la and lb with no 
magnetic diffusion (with 51 and 101 grid points in the X 
and Y axes to check the rate of numerical dissipation), Ic 
with the diffusion in three dimensions (77 = lO^^cm^/s, 51 
grid points). Id with the planar diffusion (only in the XY 
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Fig. 4. a: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model Id 
with the 2D magnetic diffusion at t — 8.5 • 10 yr (top) and at t = 9 • 10 yr (bottom) 



plane, 77 = lO^^cm^/s, 51 grid points) and le with the pla- 
nar diffusion {rj = 3.6-lO^^cm^/s, 101 grid points). The re- 
sults concerning the evolution of a magnetic field structure 
are discussed, as an example, for the case lb and Id. The 
rest of computations for the model I are presented in the 
discussion about magnetic energy density time evolution, 
later in this paragraph. In Fig. la-d the distribution of the 
stellar and gaseous disk for the model I is presented for 
the same evolutionary stages. We can see that the whole 
galactic structure is clearly visible for stars and for gas. 
The system rotates in the anti-clockwise direction. The 
evolution of magnetic field for the case lb is presented in 
Fig. 3a (top) for 8-10^ yr, where the rotation of completely 
frozen-in magnetic field, which follows the motions of gas 



is shown. The magnetic field lines are aligned with the 
bar and spiral arms, however in the inner layers of both 
arms a small reversals of magnetic field vectors (shown in 
the galactic plane) could be seen (Fig. 3a, top-left). The 
main magnetic arms are distributed slightly outside the 
main molecular arms. The magnetic intensity distribu- 
tion in the galactic plane (Fig. 3a top-right, gray plot) 
shows that the magnetic spirals, as well as spiral arms 
from molecular clouds do not form a single, well described 
feature. The magnetic arms arc dispersed into two, three 
parts, where the direction of vectors could be opposite to 
the initial field configuration. There could be also seen 
magnetic bridges between two closely distributed mag- 
netic arms (Fig. 3a, top- right). The magnetic field vectors 
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Fig. 4. b: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model Id 
with the 2D magnetic diflfusion at t = 9.5 • 10** yr 



start to be mixed at the ends of the bar (Fig. 3a, top- left). 
For 8.5 • 10^ yr (Fig. 3a, bottom) the regions, where field 
direction changes abruptly, are also present. The magnetic 
field configuration begins to be more complicated than the 
velocity field in the same moment of time, what is present 
also in the magnetic evolution performed in the 2D ex- 
periment (Otmianowska-Mazur et al. 1995). The regions 
with magnetic vector reversal are distributed in the inner 
layers of the magnetic arms (Fig. 3a, bottom- left). The 
magnetic intensity maxima are visible again in the spiral 
arms and along the bar (Fig. 3a, bottom-right) and the re- 
versal regions show significantly lower intensity. For 9 • 10^ 
yr (Fig. 3b, top) the rotation of the structure follows, the 
galaxy rotates once, field is more and more twisted and 
two additional spiral arms, going in the opposite direction 
than main spirals, form. The magnetic strength, as it is 
seen in Fig. 3b top- right (gray plot), is highest in the dy- 
namical structure, however regions with reversals are here 
better visible than for t = 8.5 • 10* yr, especially for the 
upper magnetic arm. The last graph for 9.5T0* yr presents 
regions where even three reversals of magnetic field vec- 
tors could be seen (Fig. 3b, bottom- left) . The magnetic 
configuration shows two arms in the upper part of the fig- 
ure with opposite vector direction. The intensity of the 
arm going according to the input magnetic direction is 
certainly higher, which is visible in the magnetic inten- 
sity gray plot (Fig. 3b, bottom-right). In the lower part of 
Fig. 3b, bottom-left, only one magnetic arm is visible with 
normally directed vectors; the second arm shows oppo- 
site direction. The magnetic evolution with high Reynolds 
numbers shows that the influence of a dynamical galactic 



evolution on a magnetic field history is much more compli- 
cated than we expected. The magnetic field configuration, 
as well as the velocity field, rotates in the anti-clockwise 
direction in the whole evolutionary time. The numerical 
simulations also present that the magnetic field vectors 
are mixed at both ends of a bar (e.g. Fig. 3a, bottom- left). 
The analysis of gas density maxima distributions shows 
that maxima of magnetic intensity are distributed slightly 
outside the molecular spirals. 

The presence of the magnetic diffusion changes the pic- 
ture significantly (case Id, Fig. 4). For the presentation we 
choose the case with the highest value of the coefficient 
rj — lO^^cm^/s. At the first evolutionary stages the mag- 
netic field structure is not influenced strongly by the dif- 
fusion. The situation changes for 8.5 • lO^yr (Fig. 4a, top), 
when the difference between configurations with and with- 
out dissipation becomes extremely large. Figure 4a shows 
that the magnetic structure in the galactic plane is not 
so aligned with the bar as for the model lb (Fig. 3a, bot- 
tom) . The magnetic arms are much broader and big parts 
of them are completely dispersed. The main magnetic 
arm which survives is that one with magnetic vectors di- 
rected opposite to the initial direction of the magnetic field 
(Fig. 4a, top) and is certainly not connected with dynam- 
ical structure of the model galaxy. For 9 • 10* yr the mag- 
netic intensity maxima are again more dissipated, leav- 
ing only two arms with opposite directions of the vectors 
(Fig. 4a, bottom) and with similar intensity (Fig. 4a, right) 
and rather weak magnetic arm on the left-hand side of 
the figure. Figure 4b shows magnetic intensity (right) and 
magnetic field vectors (both in the galactic plane) (left) 



K. Otmianowska-Mazur et al.: On 3D galactic magnetic field evolution 



11 



for 9.5 • 10^ yr. At this stage of evolution the magnetic en- 
ergy density is extremely low due to very strong diffusion. 
Almost the whole structure is dispersed, however in the 
upper part of the figure there could be seen two magnetic 
arms with the field reversal. The numerical experiment 
with the high value of the diffusion coefficient results in 
much more smoothed magnetic structures (Fig. 4) than for 
the simulations lb (Fig. 3). The intensity of regions, where 
magnetic vectors are distributed randomly, decreases fast 
and causes higher contrast between smoothed and random 
vector areas (Fig. 4). 



Fig. 5. The topolines of the velocity field intensity and mag- 
netic intensity as a gray plot at t = 9 ■ 10* yr for the model 
lb with no magnetic diffusion. 



Fig. 6. a: The topolines of the magnetic field component Bz 
at t = 9.0 • 10* yr in the galactic plane for the model lb. 



In order to compare the structure of the magnetic field 
with the velocity field intensity Fig. 5 presents the topo- 
lines of the velocity field intensity over the magnetic inten- 
sity gray-plots for the case without the diffusion (model 
lb, t — 9 ■ 10^ yr). The basic velocity minimum is situ- 
ated in the galactic center, then two maxima are visible 
under and above the galactic bar. The magnetic maxima 
connected with the bar coincide to some extent with the 
velocity maxima, however the magnetic spirals are dis- 
tributed slightly outward from the galactic center. The 
spiral arms are not observed in the topolines map, which 
means that velocity field gradients connected with them 
are not high. 

The magnetic configuration is examined not only for 
the galactic plane but also for levels above and under the 
XY plane. We have found that the magnetic field forms 
the configuration depending on the z coordinate. In Fig. 6a 



Fig. 6. 

plane. 



b: The same as in Fig. 6a but 350 pc above the galactic 
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Fig. 6. c: The same as in Fig. 6a but 350 pc under the galactic 
plane. 



Fig. 7. Time evolution of the magnetic energy density E 
normalised to the initial energy density Eo for the following 
models: la and lb with no magnetic diffusion for 51 and 101 grid 
points, respectively, Ic with 3D magnetic diffusion coefficient 
lO'^^ cm^/s. Id with planar diffusion coefficient 10^^ cm^/s and 



le with planar diffusion coefficient 3.6 • 



7s 



the topolincs of the B^ component for the galactic plane 
(level=ll) at time step 9 • 10^ yr is presented for the ex- 
periment lb. In Fig. 6b and c we present parallel cuts 
of these B^ component 350 pc above (level=4) and un- 
der (level=17) the galactic plane. Comparing Fig 6a with 
Fig. 6b it is visible that the magnetic field disappears in 
certain regions. Above the galactic plane the situation 
is just opposite (Fig. 6c). Regions vanishing at level=4 
now are clearly visible, but a few next regions disappear. 
The edge-on-view of the model galaxy would give different 
structures above and under the galactic plane. 

The maps of the topolines of the B^ component (Fig. 6) 
inform us also that the magnetic field lines wave in the 
Z direction going along the spiral arms and bar. The min- 
ima and maxima of the B^ component go one by one along 
these structures changing their sign from plus to minus 
(Fig. 6). The amplitude of these waves is not high - the 
value of B^ is of one order smaller than the value of i?x 
and By, but such waving is visible in most time steps of 
the experiments with and without diffusion. 

Figure 7 presents the time evolution of magnetic en- 
ergy density (E) normalised to the initial energy density 
Eq for the five experiments which have been performed 
for the model I galaxy. The figure legend is in agreement 
with the description of the model I in Table 2. In order to 
solve the problem of the numerical diffusion two experi- 
ments with no additional magnetic dissipation have been 
made. In the case la two times larger distance between 
grid points have been applied. The simulations lb result 
in significantly higher increase of magnetic energy density 
than for the case with larger step. It means that numerical 
dissipation due to truncation errors is higher for the bigger 
grid interval. We estimate that in the case la the numerical 
diffusion is of order 10^^ cm^/s and about one order less for 
the case lb. For both simulations the growth of magnetic 
energy density is present, however the increase in the ex- 
periment lb is certainly higher. The initial three maxima 
of E are followed by quick decaying of its value due to the 
effects of numerical diffusion. The magnetic energy evo- 
lution for the model lb could be compared with the case 
le with the same number of grid points and grid interval. 
The simulations le have been done for the planar coeffi- 
cient of diffusion 77 = 3.6 • 10^^ cm^/s. Figure 7 shows that 
a small increase of magnetic energy density is still present 
in this case, however after very short evolutionary time 
(8 • 10^ yr) a fast decreasing appears. The models Ic and 
Id are prepared for 51 grid points in the X and Y directions 
and could be compared with the case la. The simulations 
Ic and Id have been done for three-dimensional diffusion 
77 = 10^^ cm^/s and for the planar dissipation 10^^ cm^/s, 
respectively. For both cases the magnetic energy density 
E simply decreases with the evolutionary time due to the 
numerical and physical diffusion. The diminishing of E is 
faster for the computations Id with the biggest magnitude 
of the diffusion coefficient. The decaying of E is connected 
with magnetic energy outflow through X and Y borders of 
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Fig. 8. a: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model Ila 
with no magnetic diffusion at f — 3 • 10* yr (top) and at f = 4.5 • 10* yr (bottom) 



the model galaxy, while for the three-dimensional diffusion 
(Ic) outflow is present mainly in the Z direction. 

The maxima of the energy density evolution in case 
la and lb can be clearly identified as driven by the differ- 
ent stages of bar evolution. At timestep 76 the bar starts 
to grow, reaching its maximum strength and length at 
timestep 81 and then it shrinks slowly (up to timestep 90) 
followed by an accelerated shrinkage. 

All of our simulated models result in the decay of mag- 
netic energy density due to the numerical and/or physical 
diffusion. Each model suggests the amplification of mag- 
netic field in the arm due to compressional and shearing 
motions of gas induced by the bar and arms. The mod- 
els with the diffusion working mainly in the XY direction 
(model 1 e) show much slower decline of the energy density. 



4.4- The 3D evolution of magnetic field: model II 

The calculations of magnetic field evolution with gravita- 
tionally weak bar disturbance (model II) have been per- 
formed for earlier time stages than in the model I. 

In order to illustrate high dynamics of the simula- 
tion process, the evolution of stars and molecular clouds 
is presented in (Fig. 2a-c) for three evolutionary stages: 
t = 2 • 10^ yr, i = 3 • 10^ yr and t = 5 • 10^ yr. The 
magnetic field behaviour in this case is shown for later 
moments of time (Fig. 8). Figure 8 presents the gray plots 
of the intensity and the vectors of the magnetic field in 
the galactic plane, for the model Ila (see Table 2), with- 
out the magnetic diffusion. Figure 8a, top, corresponds to 
the evolutionary stage of t = 3 • 10® yr. Fig. 8a, bottom, 
to t = 4.5 • 10® yr, and Fig. 8b to f = 5 • 10® yr. Al- 
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Fig. 8. b: Magnetic field vectors (left) and intensity of magnetic field (gray plot, right) in the galactic plane for the model Ila 
with no magnetic diffusion at t = 5-10 yr 



ready for the first time of evolution, the initially circular 
magnetic field configuration is transformed into a struc- 
ture with a bar and spiral arms. The direction of vectors is 
mainly anti-clockwise, however in the inner region of the 
upper arm (Fig. 8a, top-left) the small area with vector re- 
versal appears. The magnetic intensity gray-plot (Fig. 8a, 
top-right) shows nonuniform distribution of the magnetic 
arms, which consist of two, three stripes due to a similar 
behaviour of the dynamical arms (Fig. 2b). According to 
dynamical features, after 1.5 • 10* yr (Fig. 8a, bottom) the 
magnetic arms become more narrow. The intensity along 
the arm (Fig. 8a, bottom-right) is still nonuniform form- 
ing some kind of bifurcation regions. The magnetic field 
waves in the XY plane in the lower magnetic arm (Fig. 8a, 
bottom-left). The region with opposite direction of mag- 
netic vectors inside the upper spiral arm is again visi- 
ble, however its intensity is rather weak (Fig. 8a, bottom- 
right). The magnetic field vector mixing is present in both 
time steps in Fig. 8a. Figure 8b illustrates that magnetic 
field structure is different from the molecular gas distri- 
bution (Fig. 2c). In the lower part of the figure the mag- 
netic bridge between two arms is present. The intensity 
map (Fig. 8b, left) presents that the magnetic arms are 
nonuniform and wave in the XY plane. 

The magnetic time evolution of energy density E nor- 
malised to the initial energy density Eq is illustrated in 
Figure 9. The simulation for the model Ila results in the 
initial growing of E. After two close maxima, E decreases 
with the evolutionary time due to numerical dissipation. 
In the case lib with 3D magnetic diffusion the magnetic 
energy density simply diminishes with time due to physi- 
cal as well as numerical diffusion. 



Calculations for model II result also in magnetic field 
mixing and waving of the magnetic lines of force in the Z 
direction. 

The magnetic field diminishing in the central part in 
all experiments is connected with assumed lower magnetic 
strength in this region and with numerical dissipation due 
to truncation errors. 



5. Conclusions and discussion 

The evolution of the large-scale galactic magnetic field 
affected by non-axisymmetric dynamical structure like 
spirals and bars has been demonstrated using 3D nu- 
merical simulations. The problem is solved with three- 
dimensional, fully dynamical velocity field evolving in 
time, resulting from N-body particle hydrodynamics, as 
the input parameter. The calculations have been per- 
formed for two dynamical experiments with strong bar 
and clearly seen spiral structure at the end of a bar. The 
basic configuration for the initial magnetic field is axisym- 
metric structure. The case of the uniform input magnetic 
field has been demonstrated in Paper I and missed as not 
interesting in the actual project. 

During the first time steps of magnetic field evolution 
the lines of force follow closely the non-axisymmetric dy- 
namical structures. After about 5 ■ 10^ yr the situation 
changes and the magnetic fields become more complicated 
than the velocity field at the given time step. The regions 
with magnetic vectors going opposite to the main mag- 
netic spirals are visible. This fact is connected with the 
molecular gas flows in the dynamical structure. Such re- 
versals could probably explain the observed reversals of 
magnetic field vectors in our Galaxy (Rand and Lyne 
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1994). Further evolution shows more and more compli- 
cated structures of magnetic lines of force and after the 
whole evolutionary time the magnetic field is still aligned 
with the bar and spirals and the regions with reversals are 
still present. In the simulation with planar diffusion they 
are still seen after 9.5 • 10^ yr. 



Fig. 9. Time evolution of the magnetic energy density E nor- 
malised to the initial energy density Eo for the model II with 
no magnetic diffusion (dashed line) and with 3D diffusion (solid 
line) 



Moreover the magnetic field arms are not single fea- 
tures. They form two or three stripes which change their 
shape during the evolution. The magnetic intensity in 
arms is not uniform creating some kind of bifurcations 
along the arms and bar. The presented simulations re- 
sulted also in the waving in the XY plane of magnetic 
field vectors (Fig. 8b, left). 

The regions with magnetic field mixing are present at 
both ends of the bar. This fact results from the rotation of 
the dynamical structure. The end of the bar is the place 
where the field aligned along two opposite sides meets (see 
Figs. 3b, bottom & 4b). The magnetic field vectors have 
opposite directions going along two bar sides and during 
rotation they meet at the bar ends. 

Both obtained facts: magnetic lines alignment along 
the spiral arms and mixing of the magnetic vectors at the 
bar ends could explain the characteristic radio features ob- 
served in the barred galaxy M83 (Neininger et al. 1991), 
however the problem needs further study. The model of 
polarised emission with the telescope beam and the Fara- 
day rotation for resulting magnetic field structures is in 



preparation. The mechanism of magnetic field tangling at 
the ends of a bar is normally explained by enhanced tur- 
bulent motions of molecular clouds (e.g. Neininger et al. 
1991). 

The non-symmetrical effect of the levels from above 
and under the galactic plane shows that the magnetised 
galactic disk is not uniform along the Z-axis. The magnetic 
field changes its configuration from the galactic plane to 
the z-edges of the disk. The magnetic field structure is also 
highly different in the part above and under the galactic 
plane. There is certainly no mirror symmetry according to 
this plane as it is normally used in a number of simulations 
(e.g. Brandenburg et al. 1992). 

The mechanism of magnetic field waving in the Z-axis 
constitutes the next interesting result from our simula- 
tions. The Bz magnetic component obtained from our cal- 
culations is of one order smaller than B^ and By . The mag- 
netic field observations show that the magnetic field in the 
direction perpendicular to the galactic plane is not high 
(Dumke et al. 1995). Such behaviour of magnetic lines of 
force has been found observationally in M31 (Beck et al., 
1989) and probably could help to initiate development of 
the Parker instability in the spiral arms and bar (Hanasz 
and Lesch 1993). 

Finally, we discuss the time-evolution of the magnetic 
energy density Wb — E/Eo. We have already mentioned 
that the decline of Wb is mainly related to the numerical 
diffusion which could be even of 10^^ cm^/s order in the 
case of big grid interval. 

Our simulations clearly indicate the relation between 
the non-axisymmetric dynamical perturbations and the 
magnetic field evolution as long as the diffusion is not too 
strong. Any local deviation from the global velocity field 
is immediately visible in the magnetic field strength and 
structure. We emphasize that in case la and lb the differ- 
ent stages of bar dynamics can be found in the evolution 
of the magnetic field energy density Wb ■ Growth of a bar 
up to its maximum strength (time step 76 to 81) and its 
subsequent shrinkage (slow up to time step 90 faster for 
time steps > 90) ) can be easily related with the max- 
ima of the energy density time profile. The same holds 
for model II. After a first maximum is reached (3.2 • 10* 
yr), which is correlated with the appearance of the 2 arm 
spirals, Wb drops a little and then increases again. The 
second maximum appears when the bar reaches his maxi- 
mum length and strength (3.9- 10* yr). After t — 4- 10* yr 
Wb decreases due to numerical diffusion connected with 
the truncation errors. 

Obviously, the magnetic field reacts sensitively to vari- 
ations in the gas velocity, triggered by dynamical pro- 
cesses. These variations can be seen from Fig. 10 a and 
b where we plotted the surface mass density E and the 
gas rotation curves for both models at time steps where 
in both models a bar has finished its growth. The dif- 
ferences in the dynamical behaviour of both models are 
obvious. In model II there is much more structure in the 
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Fig. 10. Comparison between the two bar models, 
velocity v^ of two models 



The left panel shows the surface density E, the right panel the rotation 



rotation curve and a higher central mass concentration. 
Dynamically speaking model II is more "active" . 

This enhanced activity is easy to understand, since the 
gas constant of the simulated galaxy is higher in model II. 
Dynamical disk instabilities rely on the ratio of gas mass 
to total mass, and a higher gas contents leads to a stronger 
dynamical disk activity. 
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